A WEAK GALERKIN FINITE ELEMENT METHOD FOR 
SECOND-ORDER ELLIPTIC PROBLEMS 
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Abstract. In this paper, authors shall introduce a finite element method by using a weakly 
defined gradient operator over discontinuous functions with heterogeneous properties. The use of 
weak gradients and their approximations results in a new concept called discrete weak gradients 
which is expected to play important roles in numerical methods for partial differential equations. 
This article intends to provide a general framework for operating differential operators on functions 
with heterogeneous properties. As a demonstrative example, the discrete weak gradient operator is 
employed as a building block to approximate the solution of a model second order elliptic problem, 
in which the classical gradient operator is replaced by the discrete weak gradient. The resulting 
numerical approximation is called a weak Galerkin (WG) finite element solution. It can be seen 
that the weak Galerkin method allows the use of totally discontinuous functions in the finite element 
procedure. For the second order elliptic problem, an optimal order error estimate in both a discrete 

and L-^ norms are established for the corresponding weak Galerkin finite element solutions. A 
superconvergence is also observed for the weak Galerkin approximation. 
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1. Introduction. The goal of this paper is to introduce a numerical approxi- 
mation technique for partial differential equations based on a new interpretation of 
differential operators and their approximations. To illustrate the main idea, we con- 
sider the Dirichlet problem for second-order elliptic equations which seeks an unknown 
functions u = u{x) satisfying 

(1.1) - V • (aVii) + V • (6m) -I- cu = / in 17, 

(1.2) u = g on dft, 

where $7 is a polygonal or polyhedral domain in M'* {d = 2,3), a = {aij{x))dxd 6 
[L°°(ri)]'' is a symmetric matrix-valued function, b = (6i(a;))dxi is a vector-valued 
function, and c = c{x) is a scalar function on n. Assume that the matrix a satisfies 
the following property: there exists a constant a > such that 

(1.3) aC'^C < C^aC, Ve e R'^. 

For simplicity, we shall concentrate on two-dimensional problems only (i.e., d = 2). 
An extension to higher-dimensional problems is straightforward. 

The standard weak form for ()1.1|) and (|1.2|) seeks u G H^{n) such that u = g on 
dH, and 

(1.4) {aVu,Vv)-{bu,Vv) + {cu,v) = {f,v) Vv € H^{n), 
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where ('/','(/') represents the L^-inner product of (p — (j){x) and ip = ip{x) - either 
vector-valued or scalar-valued functions. Here Vu denotes the gradient of the func- 
tion u = u{x), and V is known as the gradient operator. In the standard Galerkin 
method (e.g., see [13l[7|), the trial space H^{VL) and the test space HI{VL) in (jl.4p 
are each replaced by properly defined subspaces of finite dimensions. The resulting 
solution in the subspace/subset is called a Galerkin approximation. A key feature in 
the Galerkin method is that the approximating functions are chosen in a way that 
the gradient operator V can be successfully applied to them in the classical sense. A 
typical implication of this property in Galerkin finite element methods is that the ap- 
proximating functions (both trial and test) are continuous piecewise polynomials over 
a prescribed finite element partition for the domain, often denoted by 7^. Therefore, 
a great attention has been paid to a satisfaction of the embedded "continuity" re- 
quirement in the research of Galerkin finite element methods in existing literature till 
recent advances in the development of discontinuous Galerkin methods. But the inter- 
pretation of the gradient operator still lies in the classical sense for both "continuous" 
and "discontinuous" Galerkin finite element methods in current existing literature. 

In this paper, we will introduce a weak gradient operator defined on a space of 
functions with heterogeneous properties. The weak gradient operator will then be 
employed to discretize the problem (11.41) through the use of a discrete weak gradient 
operator as building bricks. The corresponding finite element method is called weak 
Galerkin method. Details can be found in Section |4l 

To explain weak gradients, let K be any polygonal domain with interior K'^ and 
boundary dK. A weak function on the region K refers to a vector-valued function 
V = {vo^Vh} such that vq G L'^{K) and Vb & H^{dK). The first component vq can 
be understood as the value of v in the interior of K , and the second component Vb is 
the value of v on the boundary of K . Note that Vb may not be necessarily related to 
the trace of vq on dK should a trace be defined. Denote by W{K) the space of weak 
functions associated with K] i.e., 

(1.5) W{K) = {w = {vo.Vb] : «o e L''{K), Vb G H^dK)}. 

Recall that the dual of L^{K) can be identified with itself by using the standard 
inner product as the action of linear functionals. With a similar interpretation, 
for any v G W{K), the weak gradient of v can be defined as a linear functional V^v 
in the dual space of H{div, K) whose action on each q G (div, K) is given by 

(1.6) (VdW,q):=- / VQ\7-qdK+ / Vbq-nds, 

JK JdK 

where n is the outward normal direction to dK . Observe that for any v G W{K)^ 
the right-hand side of (|1.6p defines a bounded linear functional on the normed linear 
space iJ (div, K). Thus, the weak gradient V^w is well defined. With the weak gradient 
operator being employed in (|1.4p . the trial and test functions can be allowed to 
take separate values/definitions on the interior of each element T and its boundary. 
Consequently, we are left with a greater option in applying the Galerkin to partial 
differential equations. 

Many numerical methods have been developed for the model problem (|l.ip - (|1.2p . 
The existing methods can be classified into two categories: (1) methods based on 
the primary variable u, and (2) methods based on the variable u and a flux variable 
(mixed formulation). The standard Galerkin finite element methods ([I31[71[S]) and 
various interior penalty type discontinuous Galerkin methods ([H 121 El HIl [H]) are 



3 



typical examples of the first category. The standard mixed finite elements ([301 HI 111 
[SI ini [HI Uni 121] ) and various discontinuous Galerkin methods based on both variables 
( [121 [T31 [ini [TH]) are representatives of the second category. Due to the enormous 
amount of publications available in general finite element methods, it is unrealistic to 
list all the key contributions from the computational mathematics research community 
in this article. The main intention of the above citation is to draw a connection 
between existing numerical methods with the one that is to be presented in the rest 
of the Sections. 

The weak Galerkin finite element method, as detailed in Section [H is closely 
related to the mixed finite element method (see |20l [H (H [8l [TTl [24] ) with a hybridized 
interpretation of Fraeijs de Veubeke |T71 [TH] . The hybridized formulation introduces 
a new term, known as the Lagrange multiplier, on the boundary of each element. The 
Lagrange multiplier is known to approximate the original function u = u(x) on the 
boundary of each element. The concept of weak gradients shall provide a systematic 
framework for dealing with discontinuous functions defined on elements and their 
boundaries in a near classical sense. As far as we know, the resulting weak Galerkin 
methods and their error estimates are new in many applications. 

2. Preliminaries and Notations. We use standard definitions for the Sobolev 
spaces H'^{D) and their associated inner products {■t-)s,d, norms || • ||s,u, and semi- 
norms I • \s^D for s > 0. For example, for any integer s > 0, the seminorm | • \s^d is 
given by 




with the usual notation 



a=(ai,a2), \a\ = ai + a2, d"" — d^^dl 
The Sobolev norm || • \\m,D is given by 

j=o 

The space H^{D) coincides with L^{D), for which the norm and the inner product 
are denoted by || • ||d and {■,-)d, respectively. When D — i},we shall drop the subscript 
D in the norm and inner product notation. The space iJ(div; f2) is defined as the set 
of vector-valued functions on il which, together with their divergence, are square 
integrable; i.e., 

i/(div; n) = {v: V e V ■ v G L^ifl)} . 

The norm in iJ(div; D,) is defined by 

l|v|k(div;n) = (l|vf +||V.vf)^ 

3. A Weak Gradient Operator and Its Approximation. The goal of this 
section is to introduce a weak gradient operator defined on a space of functions with 




4 



heterogeneous properties. The weak gradient operator will then be employed to dis- 
cretize partial differential equations. To this end, let K be any polygonal domain 
with interior K'^ and boundary dK. A weak function on the region K refers to a 
vector-valued function v = {vo,Vb} such that vq G L^{K) and G H^{dK). The 
first component wq can be understood as the value of v in the interior of K, and the 
second component Vb is the value of v on the boundary of K. Note that Vb may not 
be necessarily related to the trace of vq on dK should a trace be well defined. Denote 
by W{K) the space of weak functions associated with K; i.e., 

(3.1) W{K) ^{v = {v^.Vb] : vo e L^K), Vb G H^dK)}. 

Definition 3.1. The dual of L^{K) can be identified with itself by using the 
standard inner product as the action of linear functionals. With a similar inter- 
pretation, for any v G W(K), the weak gradient ofv is defined as a linear functional 
V dV in the dual space of H{d\v,K) whose action on each q G H{div,K) is given by 

(3.2) (VdW,q):=-/ vo\7-qdK+ Vbq-nds, 

Jk JdK 

where n is the outward normal direction to dK . 

Note that for any v G W{K), the right-hand side of p.2p defines a bounded linear 
functional on the normed linear space i?(div, K). Thus, the weak gradient V^w is well 
defined. Moreover, if the components of v are restrictions of a function u G H^{K) 
on K^ and dK , respectively, then we would have 



wqV • qdK + / Vbq ■ nds = — / wV • qdK + uq ■ nds = / Vm • qdK. 
K JdK Jk JdK Jk 

It follows that VdW — Vu is the classical gradient of u. 

Next, we introduce a discrete weak gradient operator by defining in a poly- 
nomial subspace of iJ(div,if). To this end, for any non- negative integer r > 0, 
denote by Pr{K) the set of polynomials on K with degree no more than r. Let 
V{K,r) C [Pr{K)]^ be a subspace of the space of vector-valued polynomials of 
degree r. A discrete weak gradient operator, denoted by Vd,r, is defined so that 
^ d,rV G V{K, r) is the unique solution of the following equation 

(3.3) / Vd,rV-qdK = - voV-qdK+ Vbq-nds, \fq&V{K,r). 

Jk ' Jk JdK 

It is not hard to see that the discrete weak gradient operator Vd,r is a Galerkin- 
type approximation of the weak gradient operator by using the polynomial space 
vlK,r). 

The classical gradient operator V — {dxi,dx2) should be applied to functions 
with certain smoothness in the design of numerical methods for partial differential 
equations. For example, in the standard Galerkin finite element method, such a 
"smoothness" often refers to continuous piecewise polynomials over a prescribed finite 
element partition. With the weak gradient operator as introduced in this section, 
derivatives can be taken for functions without any continuity across the boundary of 
each triangle. Thus, the concept of weak gradient allows the use of functions with 
heterogeneous properties in approximation. 

Analogies of weak gradient can be established for other differential operators such 
as divergence and curl operators. Details for weak divergence and weak curl operators 
and their applications in numerical methods will be given in forthcoming papers. 
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4. A Weak Galerkin Finite Element Method. The goal of this section is to 
demonstrate how discrete weak gradients be used in the design of numerical schemes 
that approximate the solution of partial differential equations. For simplicity, we take 
the second order elliptic equation as a model for discussion. With the Dirichlet 
boundary condition (jl.2p . the standard weak form seeks u G H^{fl) such that u = g 
on dVl and 

(4.1) (aVu, Vu) - Vi;) + (cm, u) = (/, u) \lv Hl{n). 

Let 7?i be a triangular partition of the domain f2 with mesh size h. Assume that 
the partition Th is shape regular so that the routine inverse inequality in the finite 
element analysis holds true (see [H]). In the general spirit of Galerkin procedure, 
we shall design a weak Galerkin method for (|4.ip by following two basic principles: 
(1) replace H^(fl) by a space of discrete weak functions defined on the finite element 
partition Th and the boundary of triangular elements; (2) replace the classical gradient 
operator by a discrete weak gradient operator T/d,r for weak functions on each triangle 
T. Details are to be presented in the rest of this section. 

For each T E Th, Denote by Pj{T^) the set of polynomials on T° with degree no 
more than j, and Pi{dT) the set of polynomials on dT with degree no more than i 
(i.e., polynomials of degree H. on each line segment of dT). A discrete weak function 

V = {wo,W6} on T refers to a weak function v — {vQ^Vb] such that G PjC^^) ^^'^ 
Vb G Pi{dT) with j >0 and£> 0. Denote this space by W{TJ,£), i.e., 

W{TJJ) := {v = {vo,Vb} : «o £ ^j(T°),Wb £ Pi{dT)} . 

The corresponding finite element space would be defined by patching W{T,j,£) over 
all the triangles T £ Th- In other words, the weak finite element space is given by 

(4.2) ShUJ)-^{v^{vo,Vb}: {vo,Vb}\T e WiT,j,e),yT e %} ■ 

Denote by S^{j,£) the subspace of Sh{j,£) with vanishing boundary values on 
i.e., 

(4.3) Slij,£) {v = {vo,Vb} G Shij,£),Vb\oTndn = 0, VT G %} . 

According to p.3p . for each v — {vo,Vb} G Sh{j,£), the discrete weak gradient of 

V on each element T is given by the following equation: 

(4.4) / Vd.rV ■qdT = - / vqV ■ qdT + / Vbq ■ nds, Vg G V{T, r). 
Jt Jt Jot 

Note that no specific examples of the approximating space V(T, r) have been men- 
tioned, except that V{T, r) is a subspace of the set of vector- valued polynomials of 
degree no more than r on T. 

For any w,v G Shij,£), we introduce the following bilinear form 

(4.5) a{w,v) = {a\7d,rW, Vd,rv) - {buo,Vd,rv) + (cuo,wo), 
where 

{aS/d,rW, Vd,rV) = / aVd.rW ■ \Jd,rVdVL, 

Jn 

(bwo, Vd.rV) ^ buo ■ Vd.rVdil, 

Jn 

{cwo,vo) = / cwovodil. 
Jn 
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Weak Galerkin Algorithm 1 . A numerical approximation for U.l\) and U.3\) 
can be obtained by seeking uj^ — {uq,ui,} G Sii{j,£) satisfying ui, — Qf,g on dfl and 
the following equation: 

(4.6) a{uh,v) = if, vo), yv = {«o,«b} € 

where Qbg is an approximation of the boundary value in the polynomial space Pi{dTr\ 
drt). For simplicity, Qbg shall be taken as the standard projection for each bound- 
ary segment; other approximations of the boundary value u = g can also be employed 
in J/^.gp . 

5. Examples of Weak Galerkin Method with Properties. Although the 
weak Galerkin scheme (j4.6|) is defined for arbitrary indices j, £, and r, the method can 
be shown to produce good numerical approximations for the solution of the original 
partial differential equation only with a certain combination of their values. For one 
thing, there are at least two prominent properties that the discrete gradient operator 
'Vd,r should possess in order for the weak Galerkin method to work well. These two 
properties are: 

PI: For any v G Sh{j, t), if 'V d,rV = on T, then one must have v = constant on 
T. In other words, vq = vi, = constant on T; 

P2: Let u G H™'{Vl){m > 1) be a smooth function on i7, and QhU be a certain 
interpolation/projection of u in the finite element space Sh{j,£)- Then, the 
discrete weak gradient of QhU should be a good approximation of Vm. 

The following are two examples of weak finite element spaces that fit well into 
the numerical scheme (|4.6|) . 

WG Example 5.1. In this example, we take i = r = and V{T,j+l) = 

[Pj+i(T)] , where j > is any non-negative integer. Denote by Sh{j,j + 1) the 
corresponding finite element space. More precisely, the finite element .space S';i(j, j + 
consists of functions v = {vi3,Vb} where vq is a polynomial of degree no more than j in 
r", and Vb is a polynomial of degree no more thanj-\-l ondT. The space V{T,r) used 
to define the discrete weak gradient operator Vd.r in Ii4-4\ l given as vector-valued 
polynomials of degree no more than j 1 on T. 

WG Example 5.2. In the second example, we take £ = j,r = j -\- 1, and 

V{T,r = i + 1) = [Pj{T)]^ + Pj{T)x, where x = {xi,X2)'^ is a column vector and 
Pj{T) is the set of homogeneous polynomials of order j in the variable x. Denote by 
Sh{j,j) the corresponding finite element space. Note that the space V{T,r) that was 
used to define a discrete weak gradient is in fact the usual Raviart- Thomas element 
\20l of order j for the vector component. 

Let us demonstrate how the two properties PI and P2 are satisfied with the 
two examples given as above. For simplicity, we shall present results only for WG 
Example 15.11 The following result addresses a satisfaction of the property PI. 

Lemma 5.1. For anyv — {woi^fc} G W{T,j,j+\), let\7 d.j+iv be the discrete weak 
gradient of v on T as defined in with V{T,r) = [Pj+i{T)]^ . Then, Vdj+iw — 

holds true on T if and only if v — constant (i.e., vq — Vb — constant). 

Proof. It is trivial to see from (14.41) that ii v — constant on T, then the right-hand 
side of (|4.4I) would be zero for any q e V{T, j + 1). Thus, we must have Vd.j+iv = 0. 
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Now assume that Vdj+iv = 0. It follows from (|4.4p that 



(5.1) 



/ VqV- qdT + [ vtq- nds = 0, Vg G V{T,j + 1). 

JT JdT 



Let vq be the average of uo over T. Using the results of [11 , there exists a vector- 
valued polynomial qi G V{T,j + 1) = [Pj+i(r)]^ such that qi ■ n = on dT and 
V • gi = Uq — vo- With g = in (jS.ip . we arrive at /^(wo — vo)^dT = 0. It follows 
that = Wo, and (|5.ip can be rewritten as 



(5.2) / {vb-VQ)q-nds^O, yqeV{T,j + 1). 

JdT 

Now since vt — vq G Pj+i(3T), then one may select a q G V{T,j + 1) = [Pj+i(T)]^ 
such that 

/ (j)q-nds^ ( <j)ivb-Vo)ds, G P^+i (<9T), 

JdT JdT 

which, together with (|5.2p and (f) = Vh — vq yields 

{vb - VQ)^ds = 0. 

dT 

The last equality implies Vb = vq = constant, which completes a proof of the lemma. 
□ 

To verify property P2, let u G H^{T) be a smooth function on T. Denote by 
QhU — {Qou, Qbu} the projection onto Pj{T^) x Pj+i{dT). In other words, on 
each element T , the function Qqu is defined as the projection of u in (T) and 
on 9T, QbU is the projection in Pj+i{dT). Furthermore, let Rh be the local 
projection onto V{T,j + 1). According to the definition of V^.j+i, the discrete weak 
gradient function ^d,j+iiQhu) is given by the following equation: 

(5.3) f Vd,j+i{Qhu)-qdT=- f iQou)V-qdT+ f {Qbu)q-nds, Vg G y(X,j+l). 

JT JT JdT 

Since Qq and Qb are L^-projection operators, then the right-hand side of (|5.3p is given 

by 

- / {Qou)V ■ qdT + / {Qbu)q ■ nds = - / uV • qdT + uq ■ nds 

JT JdT JT JdT 

= / (Vu) -qdT = I [Rh^u) ■ qdT. 

JT JT 

Thus, we have derived the following useful identity: 

(5.4) Vd.j+i{Qhu)^Rh{yu), yueH\T). 

The above identity clearly indicates that ^d,j+i{Qhu) is an excellent approximation 
of the classical gradient of u for any u G H^{T). Thus, it is reasonable to believe that 
the weak Galerkin finite element method shall provide a good numerical scheme for 
the underlying partial differential equations. 
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6. Mass Conservation of Weak Galerkin. The second order elliptic equation 
can be rewritten in a conservative form as follows: 

y ■ q + cu — f, q ^ — aVu + bu. 

Let T be any control volume. Integrating the first equation over T yields the following 
integral form of mass conservation: 



/ q ■ nds + 


I cudT — 


IdT 


Jt Jt 



(6.1) 

We claim that the numerical approximation from the weak Galerkin finite element 
method for (jl.ip retains the mass conservation property (|6.1I) with a numerical flux qh- 
To this end, for any given T G 7h, we chose in (14. 6p a test function v = {vo,Vi, — 0} 
so that Wo = 1 on T and vq — elsewhere. Using the relation (14. 5p . we arrive at 



(6.2) / aVd.rW/i • \/ d.rvdT — I buo ■ V d.rvdT + / cu^dT = I fdT. 
Jt Jt Jt Jt 

Using the definition (|4.4p for Vd,r, one has 

/ a\/d,rUh-yd,rvdT^ I Rh{aV d,rUh) ■ y d,rVdT 
Jt Jt 

= - / V • Rh{aV d.rUh)dT 
Jt 

(6.3) =-/ Rh{a\7 d,rUh) ■ nds 

JdT 



and 



/ buQ ■ Vd,rvdT = / Rh{buo) ■ V d,rvdT 
Jt Jt 

= - [ V ■Rh{buQ)dT 
Jt 



(6.4) = - / Rh{buo)-nds 

JdT 

Now substituting and ((Ol) into (lO)) yields 

(6.5) / Rh {~aVd,rUh + buo) ■ nds + / cuodT = / fdT, 
JdT Jt Jt 

which indicates that the weak Galerkin method conserves mass with a numerical flux 
given by 

qh-n = Rh {-aVd.rUh + bua) ■ n. 

The numerical flux q^ ■ n can be verified to be continuous across the edge of each 
element T through a selection of the test function v — {vq, Vf,} so that vq = and Vb 
arbitrary. 
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7. Existence and Uniqueness for Weak Galerkin Approximations. As- 
sume that w/j is a weak Galerkin approximation for tlie problem (|l.ip and (|1.2p arising 
from (|4.6p by using the finite element space Sh{j,j + 1) or Sh{j,j)- The goal of this 
section is to derive a uniqueness and existence result for Uh- For simplicity, details 
are only presented for the finite element space Sh{j, j + 1); the result can be extended 
to Sh{j,j) without any difficulty. 

First of all, let us derive the following analogy of Garding's inequality. 

Lemma 7.1. Let Sh{j,£) be the weak finite element space defined in |^.^[ ) and 
a(-, •) be the bilinear form given in |^.5[ ). There exists a constant K and ai satisfying 

(7.1) a{v,v)+K{vQ,VQ) >ai(||Vd,.wf + ||«of), 
for all V G Sh{j,i)- 

Proof. Let Bi = ||&||l=°(J2) and B2 — ||c||L=o(n) be the L°° norm of the coefficients 
b and c, respectively. Since 

\{bV0,Vd,rV)\ < Si||Vd,,.«|| \\vo\\, 

\{cvo,vo)\ < 

then it follows from (|4.5p that there exists a constant K and ai such that 

a{v,v)+K{vo,vo) > a||Vd.rvf - Si||Vd,r«|||ko|| + (K - S2)||vof 
>a,{\\Vd,rvf + \\vof), 

which completes the proof. □ 

For simplicity of notation, we shall drop the subscript r in the discrete weak 
gradient operator V^.r from now on. Readers should bear in mind that refers to 
a discrete weak gradient operator defined by using the setups of either Example 15.11 
or Example 15.21 In fact, for these two examples, one may also define a projection Ilh 
such that Il/iq e iJ(div, il), and on each T ^ Th, one has Il/iq e V{T, r = j + 1) and 
the following identity 

(V • q, vo)t = (V • n;,q, vo)t, Vwo e P,(r°). 

The following result is based on the above property of II/j . 
Lemma 7.2. For any q G H{div,n), we have 

(7.2) J2 • «o)^ = J2 (n„q, Vdv)T, 
Ten TeTh 

for all V = {vo,vi,} S Sf^{j,j + 1). 

Proof. The definition of Ilh and the definition of Vcjw imply , 

^ (-V • q, vo)t ^ (^'^ ' "o)"^ 
Ten Ten 

= ^ (n;,q,Vd7;)T. 
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Here we have used the fact that H^jq • n is continuous across each interior edge and 
vi, = on dil. This completes the proof. □ 

Lemma 7.3. For u G H'^^^{Vt) with s > 0, we have 

(7.3) \\Ilh{aVu) - aVd{Qhu)\\ < Ch'\\u\\i+s, 

(7.4) \\\7u-\7d{Qhu)\\<Ch'\\u\\i+,. 

Proof. Since from (|5.4I) we have Vd{Qhu) — i?;i(Vu), then 

||n;,(aVM) - aV4Qhu)\\ = \\Uh{aVu) ~ aRh{Vu)\\. 

Using the triangle inequahty and the definition of Ilh and Rh , we have 

||n,,(aVw) - aRhiS7u)\\ < ||n,i(aVu) - aVw|| + ||aVu - aRhiS7u)\\ 
<Ck^\u\\i+s. 



The estimate (|7.4p can be derived in a similar way. This completes a proof of the 
lemma. □ 

We are now in a position to establish a solution uniqueness and existence for the 
weak Galerkin method (|4.6p . It suffices to prove that the solution is unique. To this 
end, let e G S'^{j,j + 1) be a discrete weak function satisfying 

(7.5) a{e,v)^0, ^ {vo,Vb} € S°iij, j + 1). 

The goal is to show that e = by using a duality approach similar to what Schatz 
[25] did for the standard Galerkin finite element methods. 

Lemma 7.4. Let e = {eg, e^} e S^{j,j + 1) be a discrete weak function satisfying 
( [7.5[ ). Assume that the dual of with homogeneous Dirichlet boundary condition 

has the H^^^ regularity (s £ (0,1]^. Then, there exists a constant C such that 

(7.6) lleoll <C/i^||V,e||, 

provided that the mesh size h is sufficient small, but a fixed constant. 

Proof. Consider the following dual problem: Find w € H^{il) such that 

(7.7) - V • (aVw) - b ■ Vw + cw = eo ini} 

(7.8) w^O on dfl, 

The assumption of regularity implies that w G H^^^{fl) and there is a constant 

C such that 

(7.9) ||u.||i+. <C||eo||. 

Testing (|7.7p against cq and then using (|7.2p lead to 

||eo|P = (-V • (aVw), cq) - (6 • Vw, eo) + {cw, cq) 
= {UhiaVw), Vdc) - (Vw, 6eo) + (cw, cq) 
= (n,,(aVii;) - aVd{Qhw), V<je) + {aVdiQhw), V^e) 

~{Vw~Vd{Qhw), beo))-iVdiQhw), be^) 

+ (cw - c((3ow), eo) + (Qqw, ceo). 
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The sum of the second, forth and sixth term on the right hand side of the above 
equation equals a{e,Qhw) — due to (|7.5p . Therefore, it fohows from Lemma [7.31 
that 

lleof = (n,,(aV«;) ~ aVd{Qhw), V^e) - {Vw - Vd(Q/»u)), beo) 
+{c{w - Qow), eo) 
< Ch'\\w\\i+s (llVdell + lleoll) + Ch\\w\\i ||eo||. 

Using the i7^+*-regularity assumption (|7.9|) . we arrive at 

||eof<C/^1|eo||(||V,e|| + lleoll), 

which leads to 

lleoll <C/i^(||Vde|| + lleoll). 

Thus, when h is sufSciently small, one would obtain the desired estimate (j7.6p . This 
completes the proof. □ 

Theorem 7.5. Assume that the dual of with homogeneous Dirichlet bound- 
ary condition has H^^^ -regularity for some s £ (0,1]. The weak Gakerkin finite 
element method defined in ( [4.6] ) has a unique solution in the finite element spaces 
Shij,j + 1) and Sj{j,j) if the meshsize h is sufficiently small, but a fixed constant. 

Proof. Observe that uniqueness is equivalent to existence for the solution of (|4.6|) 
since the number of unknowns is the same as the number of equations. To prove a 
uniqueness, let u^j^^ and u'f^^ be two solutions of (|4.6p . By letting e — u\p — u']^' we 
see that (|7.5p is satisfied. Now we have from the Garding's inequality (|7.ip that 

a(e,e) + A'||eo|| > ai (|| Vde|| + ||eo||) . 

Thus, it follows from the estimate (j7.6|) of Lemma [73] that 

ai (||Vde|| + lleoll) <Ci^/^'^||Vde|| 

for h being sufficiently small. Now chose h small enough so that CKh^ < Thus, 

llVdell + lleoll =0, 

which, together with Lemma |5.1[ implies that e is a constant and eo — 0. This shows 
that e = and consequently, uj^^"* = □ 

8. Error Analysis. The goal of this section is to derive some error estimate for 
the weak Galerkin finite element method (|4.6p . We shall follow the usual approach 
in the error analysis: (1) investigating the difference between the weak finite element 
approximation Uh with a certain interpolation/projection of the exact solution through 
an error equation, (2) using a duality argument to analyze the error in the norm. 

Let us begin with the derivation of an error equation for the weak Galerkin ap- 
proximation Uh and the projection of the exact solution u in the weak finite element 
space Sh{j,j + 1). Recall that the projection is denoted by QhU = {Qou,Qbu}, 
where Qo denotes the local projection onto Pj{T) and Qb is the local projection 
onto Pj+i{dT) on each triangular element T e 7ft,. Let v = {uo, Vb} £ S^Ht] + 1) be 
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any test function. By testing (ll.ip against the first component Vq and using (|7.2p we 
arrive at 

{f,vo) ^ ^ (-V • (aVu), ^;o)t + (V • [hu], wo) + (cm, wq) 
= (n,,(aVu), Vdw) - (n,,(H, Vdw) + (cu, wo). 

Adding and subtracting the term a{QiiU,v) = (aVdiQhu), V^w) ~ {b{QQu), V^u) + 
{c{Qou), Vq) on the right hand side of the above equation and then using (j5.4p we 
obtain 

(8.1) (/,wo) = {aVdiQhu), Vdv) - {bQou, V dv) + [cQqu, vq) 

+(nft(aVu) - aRhiVu), \7dv) 

-{Uhibu) - bQou, Vdv) + (c(u - Qqu), vq), 

which can be rewritten as 

a{uh,v) = a{QhU, v) + {Ilh{aVu) - aRh{Vu), V dv) 

-(n,i(6u) - bQou, Vdv) + {c{u - Qqu), vq). 

It foUows that 

a{uh - QhU, v) = (UhiaVu) ~ aRh{Vu), V dv) 

(8.2) -{Ilh{bu) - bQou, Vdv) + {c{u - Qm), wo)- 

The equation (|8.2[) shah be called the error equation for the weak Galerkin finite 
element method (|4.6p . 

8.1. An estimate in a discrete i/^-norm. We begin with the following lemma 
which provides an estimate for the difference between the weak Galerkin approxima- 
tion Uh and the projection of the exact solution of the original problem. 

Lemma 8.1. Let u e i/^(il) be the solution of U.l]) and S1.2\} . Let Uh G Sh{j,j + 
1) be the weak Galerkin approximation of u arising from {4^.&^ . Denote by en ■= 
Uh — QhU the difference between the weak Galerkin approximation and the L^ projection 
of the exaction solution u — u{xi,X2)- Then there exists a constant C such that 

Y(l|V,(e^)f + ||e^,of ) < C (||n„(aVw) - aR,,{Vu)f + \\c{u - Q^uW 

(8.3) + \\Iih{bu) - 6Qouf ) + K\\uo - Qouf. 

Proof. Substituting v in (|8.2p by Ch ■— Uh — QhU and using the usual Cauchy- 
Schwarz inequality we arrive at 

a(e?i, Ch) = (n/i(aVu) - ai?/i(Vu), Vd{uh ~ Qhu)) 

-[Xlh{bu) - bQou, Vd[uh - Qhu)) + {c{u - Qqu), uq - Qqu) 
< ||n,,(aVu) - aRhiVu)\\ \\Vdiuh - Qhu)\\ 

+ \\Uh{bu) ~ bQou\\ \\Vd(uh - Qhu)\\ + \\c{u - Qou)\\ \\uo - Qqu\\. 
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Next, we use the Garding's inequality (17.11) to obtain 

ai(||Vd(e,Of + ||e,,,of ) < \\TlhiaVu) - ai?,,(Vu)|| ||VdK - Qhu)\\ 

+ \\Uhibu) ~ hQ^u\\ W^aiuh - Qhu)\\ 
+ ||c(w - Qqu)\\ ||mo - Qoull + K\uq - Qo^ll^ 

< ydlVdK - Qhu)f + ho - Qo"f ) 

+C (lln^(aVM) - ai?^(Vu)f + ||n,,(6M) ~ hQ^uf 
+ llc(u-gou)f) +Aluo-Qowf , 

which implies the desired estimate (|8.3p . □ 

8.2. An estimate in We use the standard duality argument to derive 

an estimate for the error Uh — QhU in the standard L? norm over domain f2. 

Lemma 8.2. Assume that the dual of the problem U.l\) and \1.2\) has the 
regularity. Let u £ H^(fl) be the solution U.l\) and and Uh be a weak Galerkin 

approximation of u arising from \4.b^ by using either the weak finite element space 
Sh{j,j + 1) or Sh{j,j)- Let QhU be the L^ projection of u in the corresponding finite 
element space (recall that it is locally defined). Then, there exists a constant C such 
that 

\\Qqu~u4< Ch' {h\\f-Qof\\ + \\Vu~RhiS7u)\\ + \\aVu - Rhia\7u)\\ + \\u-Qou\\ 
+ \\bu-Rh{bu)\\ + ||cw-goM|| + \\VdiQhU-Uh)\\), 

provided that the meshsize h is sufficiently small. 

Proof. Consider the dual problem of (|l.ip and p.2p which seeks w £ ^^o(^^) 
satisfying 

(8.4) — V • (aVw) — b ■ Vw + cw = Qqu — uq in 51 

The assumed H^^^ regularity for the dual problem implies the existence of a constant 
C such that 

(8.5) \\w\\i+s < CWQ^u - u4. 
Testing (j8.4p against Qqu ~ uq element by element gives 

IQo" - wolP = (-V • (aVw), Qou - uq) - {b ■ S7w, Qqu - uq) + {cw, Qqu - uq) 

(8.6) =1 + 11 + III, 

where /, //, and /// are defined to represent corresponding terms. Let us estimate 
each of these terms one by one. 

For the term /, we use the identity (|7.2p to obtain 

I = (-V ■ {aVw),Qou - Mo) = (n/,.(aVw), Vd(Q/iW - u/i)). 

Recall that VdiQhu) = Rh{\7u) with Rh being a local L^ projection. Thus, 

/ = {Uh{a\7w), VdiQhU - Uh)) = {UhiaVw),RhVu - VdU,,) 
= (UhiaVw), Vu - VdUh) 

(8.7) = (n,i(aVw) - aVw, Vu - V dUh) + {aVw, Vu - VdUh)- 



14 



The second term in the above equation above can be handled as follows. Adding and 
subtracting two terms (oVdQhW, VdUh) and {a{Vw — Rh'Vw), Vu) and using the fact 
that VdiQhu) — Rh{^u) and the definition of Rh, we arrive at 

(aVw, Vu — VdUh) = (aVw, Vu) — {aVw, VdUh) 

= {aVw, Vu) - {aVdQhW, S/dUh) - (a(Vw - RhVw), V dUh) 
= {aS7w, Vu) - {aVdQhW, VdUh) ~ {a{Vw - RhVw), VdUh - Vu) 
-{a{Vw - i?/iVu;), Vu) 
(8.8) = (aVw, Vu) - {aVdQhW, VaUh) - {a{Vw - RhVw), VdU,, - Vu) 

-(Vw - RhVw, aVu ~ Rh{aVu)). 

Substituting (H^) into yields 

(8.9) / = {Iih{aVw) - aVw, Vu - V dUh) - (a(Vw - RhVw), V dUh - Vu) 

— (Vw — RfiVw, aVu — i?/i(aVu)) + (aVw, Vu) — {aV dQhW, V dUu)- 

For the term //, we add and subtract {y d{Qhw), b{Qou — uq)) from // to obtain 

II = —{b ■ Vw, Qqu — uq) 

= -(Vw - VdiQhw), b{Qou - Uq)) - (VdiQhw), b{Qau - uq)) 

= -(Vw - VdiQhw), b{Qou - Uq)) - {VdiQhw), bQou) + (VdiGhw), buo). 

In the following, we will deal with the second term on the right hand side of the above 
equation. To this end, we use ()5.4p and the definition of Rh and Qo to obtain 

{VdiQhw), bQou) = {VdiQhw) - Vw, bQou) + (Vw, bQou) 

= (VdiQhw) - Vw, bQau - bu) + {Vd{Qhw) - Vw, bu) 

+ {Vw, bQou ~ bu) + {Vw, bu) 
= (VdiQhw) - Vw, bQou - bu) + (Rh{Vw) - Vw, bu - Rh{bu)) 

+ {b ■ Vw - Qo{b ■ Vw), Qou - u) + (Vw, bu). 

Combining the last two equations above, we arrive at 

(8.10) // = -(Vw - VdiQhw), b{Qou - uo)) - (VdiQhw) ~ Vw, bQ^u - bu) 

-{Rh{Vw) - Vw, bu - Rh{bu)) ~ {b-Vw~ Qo{b ■ Vw), Qqu - u) 
-(Vw, bu) + {Vd{Qhw), buo). 

As to the term ///, by adding and subtracting some terms and using the fact 
that Qo is a local projection, we easily obtain the following 

/// = (cw, Qou - Uo) = (cw - cQow, Qou - uq) + [cQow, Qqu - uo) 
= (cw - cQow, Qou - Uo) + (cQow, Qou) - {cQow, uq) 
= {cw — cQow, Qou — Uo) + {cQow - cw, Qou) + (cw, Qqu - u) 

+(cw, u) - {cQow, Uo) 
= (cw — cQow, Qou — Uo) + {Qqw - w, cQqu - cu) + {Qow — w, cu - Qo{cu)) 

+{cw - Qo{cw), Qou - u) + (cw, u) - (cQo^, uo). 
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Note that the sum of the last two terms in / (see (|8.9p ). // (see (|8.10p ). and /// 
(see the last equation above) gives 

{aVw,\7u)-{a\7dQhW,VdUh) - (Vw, bu) + (VdiQhw), buo) + {cw, u) - {cQqw, uq) 
= a{u, w) - a{uh, Qhw) 
= (/, w) - (/, Qow) 
= if - Qof, w - Qow). 

Thus, the sum of /, //, and /// can be written as follows: 

WQou ~ wolP = (/ - Qof, w - Qow) + (nft(aVw) - aVw, Vu ~ V dUh) 

- {a{Vw - RhVw),\7dUh - Vu) - ((Vw - RhVw),a\7u - RhiaVu)) 

- {Vw - VdiQhw), biQou - Mo)) - (VdiQhw) - Vw, bQou ~ bu) 

- (RhiVw) - Vw, bu - Rhibu)) - {b-Vw- Qoib ■ Vw), Qqu - u) 
+ (cw — cQqw, Qqu — uo) + {Qqw — w, cQqu — cu) 

(8.11) + (Qqw - w, cu - Qo{cu)) + {cw ~ Qo{cw), Qqu - u). 

Using the triangle inequality, (|5.4p and (|8.5I) . we can bound the second term on the 
right hand side in the above equation by 

\{Uh{a\7w) - aVw, Vu - VdUh)\ < |(n,,(aVw;) - aVw, Vu - VdQhu)\ 

+ |(n,i(aVu;) - aVw, VdQhU - VdUh)\ 
< Ch' {\\Vu - RhiVu)\\ + WVdiQhU - uh)\\) WQou ~ 

The other terms on the right hand side of (jS.lip can be estimated in a similar fashion, 
for which we state the results as follows: 

\{a{Vw - RhVw),VdUh - Vu)\ < Ch' {\\Vu - Rh{Vu)\\ + \\Vd{QhU - Uh)\\) \\Qou 
\iiVw ~ RhVw),aVu ~ Rh{aVu))\ < Ch" \aVu - RhiaVu)\\ \\Qou - mo||, 
\{Vw ~ VdiQhw), b{Qnu - uo))| < Ch" \Qau - wof , 
liVdiQhw) - Vw, bQou -bu)\ < Ch" \u - Qau\\ \\Qou - wo||, 
\iRhiVw)-Vw, bu - Rh{bu))\ < Ch' \bu - Rhibu)\\ \\Qou - uqI 
1(6 • Vw - Qoib ■ Vw), Qou - u)\ < Ch" \u ~ Qqu\ \Qqu - ito||, 
\[cw - cQqw, Qqu - uo)| < Ch\\Qou - uqH^, 
\{Qow ~ w, cQou - cu)\ < Ch\\u - Qqu\\ WQqu - uqII, 
\{Qow - w, cu~ Qo(cm))| < Ch\\cu - Qo(cu)|| HQo" - "oil, 
\{cw - Qo{cw), Qqu ~ u)\ < Ch\\u ~ Qou\\ \\Qou - uo||. 

Substituting the above estimates into (IS.lip yields 

WQou-uof < Ch' {h\\f - Qo/ll + ||V« - Rh{yu)\\ + llaVu - Rh{aVu)\\ + \\u - Qou\\ 
+ \\bu - Rh{bu)\\ + \\cu - Qo{cu)\\ + \\Vd{QhU - Uh)\\ + \\Qou ~ uq\\) \\Qqu - uqII- 

For sufficiently small meshsize /i, we have 

IIQo" - ^^oll < Ch" {h\\f - QofW + \\Vu - Rh{Vu)\\ + llaVu - Rh{aVu)\\ + \\u - Qo"|| 
+ \\bu-Rh{bu)\\ + \\cu-Qa{cu)\\ + \\V d{QhU - Uh)\\) , 

which completes the proof. □ 
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8.3. Error estimates in and L^. With the results estabhshed in Lemma 
18.11 and Lemma 18. 2[ we are ready to derive an error estimate for the weak Galerkin 
approximation Uh- To this end, we may substitute the result of Lemma 18.21 into the 
estimate shown in Lemma |8. II If so, for sufficiently small meshsize h, we would obtain 
the following estimate: 

\\yd{uh - Qhu)f + \\uo - Qouf < C {\\Iih{aVu) - ai?,,(Vu)f + \\c{u - Qqu)\\'' 

+Ch^' {h^Wf - Oo/lP + ||V« - i?„(V«)f 
+ llaVu - Rh{aWu)f + \\u - Qauf 
+ \\bu - Rh{bu)f + \\cu - Qo{cu)\\^) . 

A further use of the interpolation error estimate leads to the following error estimate 
in a discrete norm. 

Theorem 8.3. In addition to the assumption of Lemma \8.2[ assume that the 
exact solution u is sufficiently smooth such that u G iJ™+^(f2) with < m < j + 1. 
Then, there exists a constant C such that 

(8.12) llVrfK - Qhu)\\ + ho - Qqu\\ < C(/i"||u|U+i + h^+'\\f - Qo/ll). 



Now substituting the error estimate (I8.12p into the estimate of Lemma 18. 2[ and 
then using the standard interpolation error estimate we obtain 

\\UH - Qhu\\ < C [h^+'Wf - Qo.f II + /i"+1k|U+i + /i^(/i"||u|U+i + h^+'\\f - Qo/ID) 
< C (/ii+ll/ - Q0./II + /^"+1|w||™+i) • 

The result can then be summarized as follows. 

Theorem 8.4. Under the assumption of Theorem \8.SX there exists a constant C 
such that 

\\uh - Qhu\\ < C {h^+'\\f - Qo/ll + /j™+1|w|U+i) , s e (0, 1], m e (0, j + 1], 

provided that the mesh-size h is sufficiently small. 

If the exact solution u of (|l.ip and (|1.2p has the H^^^ regularity, then we have 
from Theorem 18.41 that 

\\uh~Qhu\\ <c{h^+'y\\f\\j + y+'+^\\u\\,+2) 
<ch^+^+H\\f\Vj + \\u\Vj+2) 

for some < s < 1, where s is a regularity index for the dual of (|l.ll) and (|1.2p . In 
the case that the dual has a full (i.e., s = 1) regularity, one would arrive at 

(8.13) \\uh - Qhu\\ < Cy+^ (11/11, + ||u||,+2) . 

Recall that on each triangular element T*^, the finite element functions are of polyno- 
mials of order j > 0. Thus, the error estimate (|8.13|) in fact reveals a superconvergence 
for the weak Galerkin finite element approximation arising from ()4.6|) . 
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